R Notebook

rr library(Seurat) library(monocle)

Loading required package: Biobase
Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: ‘BiocGenerics’

The following objects are masked from ‘package:parallel’:

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from ‘package:Matrix’:

    colMeans, colSums, rowMeans, rowSums, which

The following objects are masked from ‘package:stats’:

    IQR, mad, sd, var, xtabs

The following objects are masked from ‘package:base’:

    anyDuplicated, append, as.data.frame, basename, cbind, colMeans,
    colnames, colSums, dirname, do.call, duplicated, eval, evalq,
    Filter, Find, get, grep, grepl, intersect, is.unsorted, lapply,
    lengths, Map, mapply, match, mget, order, paste, pmax, pmax.int,
    pmin, pmin.int, Position, rank, rbind, Reduce, rowMeans, rownames,
    rowSums, sapply, setdiff, sort, table, tapply, union, unique,
    unsplit, which, which.max, which.min

Welcome to Bioconductor

    Vignettes contain introductory material; view with
    'browseVignettes()'. To cite Bioconductor, see
    'citation(\Biobase\)', and for packages 'citation(\pkgname\)'.

Loading required package: VGAM
Loading required package: stats4
Loading required package: splines
Loading required package: DDRTree
Loading required package: irlba

rr library(dplyr)

Attaching package: ‘dplyr’

The following object is masked from ‘package:Biobase’:

    combine

The following objects are masked from ‘package:BiocGenerics’:

    combine, intersect, setdiff, union

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

rr library(ggplot2) library(ggrepel) library(GOplot)

Loading required package: ggdendro
Loading required package: gridExtra

Attaching package: ‘gridExtra’

The following object is masked from ‘package:dplyr’:

    combine

The following object is masked from ‘package:Biobase’:

    combine

The following object is masked from ‘package:BiocGenerics’:

    combine

Loading required package: RColorBrewer

rr seurobj <- readRDS(‘../output/seurat_objects/180831/10x-180831’) cds <- readRDS(‘output/monocle/180831/10x-180831-monocle-monocle_genelist_T1T2T3_T4T5_res.1.5’)

cannot open compressed file 'output/monocle/180831/10x-180831-monocle-monocle_genelist_T1T2T3_T4T5_res.1.5', probable reason 'No such file or directory'Error in gzfile(file, \rb\) : cannot open the connection

Figure 1

rr

save_plot(\../figures/figures_paper/main_figures/figure1/10x-180831_tsne_depots.pdf\, p1, base_width=6, base_height=4.5)
save_plot(\../figures/figures_paper/main_figures/figure1/10x-180831_tsne_timepoints.pdf\, p2, base_width=4, base_height=3)

save_plot(\../figures/figures_paper/main_figures/figure1/10x-180831_tsne_t1_aligned.pdf\, p3, base_width=6, base_height=4.5)

Figure 2

Figure 3

Heatmap

Heatmap

rr

load('../output/monocle/180831/heatmaps')
heatmap <- heatmaps$heatmap_logFC0.3_ncluster6_filtered_genes
gridExtra::grid.arrange(heatmap$ph_res$gtable)

rr

save_plot(\../figures/figures_paper/main_figures/figure3/BEAM_heatmap_logFC0.3_ncluster6_filtered_genes.pdf\, heatmap$ph_res$gtable, base_width=7, base_height=8)

Pseudotime plots

Pseudotime plots

rr

#new_cds <- buildBranchCellDataSet(cds, progenitor_method='duplicate', branch_states=c(2,3), branch_point=1, branch_labels=c('Metabolic', 'ECM'), stretch=T)
#saveRDS(new_cds, file='../output/monocle/180831/10x-180831-monocle-monocle_genelist_T1T2T3_T4T5_res.1.5_duplicated-progenitor_stretched-pseudotime')
#new_cds <- readRDS('../output/monocle/180831/10x-180831-monocle-monocle_genelist_T1T2T3_T4T5_res.1.5_duplicated-progenitor_stretched-pseudotime')
#new cds contains duplicated progenitor cells (one assigned to each branch) and stretched Pseudotime
df_pdata <- pData(new_cds)[c('original_cell_id', 'Branch', 'State.labels', 'Pseudotime')]

Error in pData(new_cds) : object 'new_cds' not found

rr

width=6
height=5
save_plot(\../figures/figures_paper/main_figures/figure3/pseudotime_ADIPOQ.pdf\, adipoq_pt, base_width=width, base_height=height)
save_plot(\../figures/figures_paper/main_figures/figure3/pseudotime_UCP2.pdf\, ucp2_pt, base_width=width, base_height=height)

save_plot(\../figures/figures_paper/main_figures/figure3/pseudotime_DCN.pdf\, dcn_pt, base_width=width, base_height=height)
save_plot(\../figures/figures_paper/main_figures/figure3/pseudotime_APOD.pdf\, apod_pt, base_width=width, base_height=height)

Featureplots

rr genes <- c(‘ADIPOQ’, ‘UCP2’, ‘DCN’, ‘APOD’) featureplots_noleg <- FeaturePlot(seurobj, features.plot=genes, cols.use=c(‘gray’, ‘blue’), no.legend=T, do.return=T, no.axes=T, pt.size=1, max.cutoff=c(5))

rr featureplots_leg <- FeaturePlot(seurobj, features.plot=genes, cols.use=c(‘gray’, ‘blue’), no.legend=F, do.return=T, no.axes=T, pt.size=1, max.cutoff=c(5))

rr adipoq_fp <- featureplots_noleg\(ADIPOQ + theme(plot.title=element_blank()) ucp2_fp <- featureplots\)UCP2 + scale_color_gradient(name=‘Expression’, low=‘gray’, high=‘blue’, guide=‘colorbar’, limits=c(0,5)) + theme(plot.title=element_blank(), legend.title=element_text(size=20), legend.text=element_text(size=20), legend.key.height = unit(1.3, ‘cm’))

rr save_plot(../figures/figures_paper/main_figures/Figure3_beam/featureplots_metabolic.pdf, metabolic_fp, base_width=width, base_height=height)

Error in save_plot(\../figures/figures_paper/main_figures/Figure3_beam/featureplots_metabolic.pdf\,  : 
  object 'height' not found

GO bubble plot

rr name <- ‘heatmap_logFC0.3_ncluster6’ folder <- paste(‘../tables/tables_paper/all_tables/BEAM/’, name, sep=‘’) #clusters <- c(1,2,3,4,6) clusters <- c(2,6,3) gsea <- list() genelist <- list() for (i in clusters){ gsea_cluster <- read.table(paste(folder,’/REVIGO_cleaned/‘, ’REVIGO_cluster’, i, ’_cleaned.tsv’, sep=‘’), header=T) gsea_cluster[’cluster’] <- i genes_cluster <- read.table(paste(folder,’/genelist_cluster’, i, ‘.tsv’, sep=’’), header=T) gsea[[i]] <- gsea_cluster genelist[[i]] <- genes_cluster } gsea <- do.call(, gsea) genelist <- do.call(, genelist) create_terms <- function(gsea){ return(data.frame(category=gsea\(gsea.domain, ID=gsea\)revigo.term_ID, term=gsea\(revigo.description, adj_pval=gsea\)gsea.p.value, genes=gsea\(gsea.intersection, cluster=gsea\)cluster)) } create_genes <- function(genelist){ return(data.frame(ID=genelist\(gene_short_name, logFC=genelist\)avgLogFC_State2_State3)) } terms <- create_terms(gsea) genes <- create_genes(genelist) circ <- circle_dat(terms, genes)

rr

save_plot('../figures/figures_paper/main_figures/Figure3_beam/180831_beamGO_all-terms.pd_labels3_zscore4.5.pdf', p, base_width=20, base_height=7)

Percent.mito

rr

percent_mito_fp <- FeaturePlot(seurobj, features.plot=c('percent.mito.100'), cols.use=c('gray', 'blue'), no.legend=F, do.return=T, no.axes=T, max.cutoff = c(10))$percent.mito + ggtitle('') + scale_color_gradient(name='% mitochondrial\ngene expression', low='gray', high='blue', guide='colorbar') + theme(legend.key.height = unit(1.6, 'cm'))

rr

save_plot(\../figures/figures_paper/main_figures/figure3/featureplot_PERCENT_MITO.pdf\, percent_mito_fp, base_width=width, base_height=height)

BATLAS

rr

batlas_results <- read.table('../tables/tables_paper/supplementary_tables/BATLAS/BATLAS_branch_high_res.txt', sep='\t', header=T)
batlas_results['pseudotime_decile'] <- seurobj@meta.data$pseudotime_decile[match(batlas_results$sample, seurobj@meta.data$branch_high_res)]
batlas_results['branch'] <- seurobj@meta.data$branch[match(batlas_results$sample, seurobj@meta.data$branch_high_res)]
batlas_results <- batlas_results %>% mutate(block_branch = if_else(branch == 'Progenitor', 'Beginning', 'End'))
batlas_results['brown_centered'] <- batlas_results$brown - 0.5
batlas_results['ratio'] <- batlas_results$brown/batlas_results$white

rr

p <- ggplot(batlas_results, aes(x=as.factor(pseudotime_decile), y=brown_centered, fill=branch)) +
  geom_col(position = position_dodge2(width=0.9, preserve='single')) +
  labs(x=\Pseudotime\, fill=\Branch\, y=NULL) +
  scale_fill_manual(values=colormap.branches, labels=c('L', 'U', 'P')) +
  facet_grid(cols=vars(block_branch), switch='x') + 
  scale_y_continuous(limits=c(-0.5, 0.5), breaks=seq(-0.5, 0.5, by=0.5), labels=c('100% white', '50/50', '100% brown')) +
  theme(strip.background=element_blank(), axis.text.x=element_blank(), strip.placement='outside')
p

rr

save_plot('../figures/figures_paper/main_figures/figure3/BATLAS_results.pdf', p, base_height=5, base_width=8)

SOM Figures

Clusters

TF analysis

rr

fig <- plot_grid(
  TSNEPlot(seurobj, group.by='res.0.5', pt.size=pt.size.tsne, do.label=T),
  TSNEPlot(seurobj, group.by='res.0.7', pt.size=pt.size.tsne, do.label=T),
  TSNEPlot(seurobj, group.by='res.1', pt.size=pt.size.tsne, do.label=T),
  TSNEPlot(seurobj, group.by='res.1.5', pt.size=pt.size.tsne, do.label=T),
  ncol=2, labels='auto'
)

Metabolic and ECM genes

rr

FeaturePlot(seurobj, features.plot='Pseudotime')